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Energy transfer and constrained 
simulations in isotropic turbulence 

By J. Jimenez 1 


1. Motivation 

The defining characteristic of turbulent flows is their ability to dissipate energy 
even in the limit of zero viscosity. The Euler equations, if constrained in such a 
way that the velocity derivatives remain bounded, conserve energy. But when they 
arise as the limit of the Navier-Stokes (NS) equations, when the Reynolds number 
goes to infinity, there is persuasive empirical evidence that the gradients become 
singular as just the right function of Re for the dissipation to remain non-zero and 
to approach a well defined limit (Sreenivasan 1984). It is generally believed that 

is muting value of the dissipation is a property of the Euler equations themselves, 
independent of the particular dissipative mechanism involved, and that it can be 
normalized with the large scale properties of the turbulent flow (e.g. the kinetic 
energy per unit volume u /2 /2, and the integral scale L) without reference to the 
Reynolds number or to other dissipative quantities (Batchelor 1953, §6). This is 
usually taken to imply that the low wave number end of the energy spectrum far 
from the dissipative range, is also independent of the particular mechanism chosen 
to dispose of the energy transfer. 

While these considerations have proved adequate in analyzing the solutions of 
the NS equations at different Reynolds numbers, they are recently being used in a 
more general context to predict the effect of different subgrid models in Large Eddy 
Simulations (LES) of turbulent flows. The LES equations are neither the Euler 
nor the NS equations. Because of the intrinsically finite capacity of computers, the 
Euler equations are first truncated to a relatively low number of degrees of freedom, 
generally much lower than that needed to represent the dissipative scales, and the 
missing dissipation is provided by the addition of a “subgrid” model whose goal is to 
mimic the effect on the large scales of the degrees of freedom that have been filtered 
in the original truncation. In most cases, the resulting model is very different from 
the quadratic viscosity characteristic of the NS equations, and the hope that it will 
approximate the behavior of the original high Reynolds number flow hinges on the 
belief that the dissipation independence observed empirically for the energy transfer 
mechanism of the NS equations will persist in the new systems. 

This is, therefore, a good time to re-evaluate the original assumptions and the 
units of their applicability. This will also give us a new tool to investigate the 
mechanism of the energy cascade itself and its relation to normal viscosity. It has 
always been an open question how much of NS turbulence is unique to the solutions 
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of the NS equations, and what parts of it would be shared by other nonlinear 
systems of partial differential equations. The arguments outlined above suggest that 
the quadratic nonlinearity might be essential but that the dissipative model is not. 
While physical experiments with non-Newtonian fluids are difficult and relatively 
rare numerical experiment with different viscosity models are commonplace and 
have been used for a long time to achieve higher Reynolds numbers that those 
obtainable with Direct Numerical Simulation, especially in two dimensional flows 
As argued above, most LES computations also fall into this category. Althoug 
many of these simulations have been quite successful in approximating the resu s 
of experiments, relatively few detailed comparisons are available, and the presence 
in most cases of adjustable parameters makes the agreement less compelling. 

Moreover, there is some experimental evidence that this dissipation independence 
is not complete in all cases. The addition of small concentrations of polymers in wall 
bounded flows is known to produce a dramatic decrease of skin friction (Berman 
1978) Since skin friction for a given flow rate is directly proportional to energy 
dissipation, and since polymers are thought to act oulyon small st ™ c * ures ’ 
that evidence is troubling. The same is true of nblets (Walsh 1990), which so 
modify friction even if they are small scale features. Finally, the friction coefficient 
of smooth-walled channels does not asymptote to a fixed value, even at the hig es 
Reynolds numbers observed (Dean 1978). While wall bounded flows may be different 
from isotropic or free shear flows, most conspicuously because in the wall region the 
“local” Reynolds number can always be argued to be low, these observations should 
be taken into account when extending the results of the latter into the behavior 
of the former. In this work, we will only concern ourselves with the behavior of 
isotropic periodic numerical simulations of the Euler or LES equations and wit 
their relation to NS turbulence. 

Even in this case, we have to qualify our support for the hypothesis of complete 
independence of the large scale motions from the details of the small scales. It is well 
known that a consequence of an inviscid energy cascade is the famous Kolmogorov 
(1941) E(k) ~ *r 5/3 spectral law. Since this result is independent of viscosity, 
the arguments above would make it a property of the Euler equations. But it has 
been known for a long time that those equations have at least another power law 
equilibrium spectrum E(k)~ k 2 , corresponding to energy equipartition among the 
spectral modes (Lee 1952). In physical systems, flow fields with such a spectrum are 
highly singular and are not expected to be observable since any residual viscosity, 
however small, would damp the high wave numbers and force a decreasing spectral 
tail These solutions are, however, consistent with the inviscid equations, and t ey 
reappear as soon as the possibility of a singularity is removed. In fact, as soon as 
the Euler equations are truncated to a finite number of degrees of freedom, ey 
tend spontaneously to the equipartition solution (Fig. 1 ). This is relevant to our 
discussion because it shows that the effect of dissipation in the NS equations is not 
only to fix the length scale of the viscous limit (the Kolmogorov scale 17), but to select 
the “desired” fc -5 / 3 spectrum away from the “contamination” of the equipartition 
component k 2 . It is not immediately clear whether other dissipation mechanisms 
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FIGURE 1. Evolution of the three-dimensional energy spectrum of the Euler 
equations in a 64 3 numerical box. Time increases with energy at high wave numbers. 
The initial flow field is fully developed turbulence from a NS simulation (solid line). 

are able to do the same. 

In the following sections, we present some numerical experiments on the effect 
of substituting different dissipation models into the truncated Euler equations. We 
will see that the effect is mainly felt in the “near dissipation” range of the energy 
spectrum, but that this range can be quite wide in some cases, contaminating a 
substantial range of wave numbers. In the process, we will develop a “practical” 
approximation to the subgrid energy transfer in isotropic turbulence, and we will 
gain insight into the structure of the nonlinear interactions among turbulent scales 
of comparable size, and into the nature of energy backscatter. Some considerations 
on future research directions are offered at the end. 

2. Inertial simulations 

2.1 The truncated Euler equations 

This section of the paper will deal with attempts to model the inertial spectral 
range of isotropic turbulence at infinite Reynolds number. 

Our first numerical experiment simulates spatially periodic turbulence with no 
viscosity. The numerical method is fully spectral, using primitive variables u, p, 
with dealiasing achieved by spherical masking in Fourier space and phase shifting 
(Rogallo, 1981, see description in Canuto et a/., 1987). The resolutions quoted for 
the different simulations reflect the number of real Fourier modes in each direc- 
tion before dealiasing; the number of actual useful modes in the power spectra, 
&max? is slightly less than half that number. The Fourier expansion functions are 
exp (ikjXj), |^| = 0, so that the length of the box side is always 2 ti\ The 

time stepping procedure is second-order Runge-Kutta for the nonlinear terms and 
an analytic integrating factor for the viscous ones, when present. The time step is 
automatically controlled to satisfy the numerical stability condition. Unless stated 
otherwise, all experiments are forced to achieve a statistically stationary flow, and 




174 


J. Jimenez 


the results presented are long term averages, usually taken over periods of the order 
of one or two large scale turnover times after stationarity is achieved. Forcing is 
achieved by multiplying all the modes in a spherical spectral shell by an appropriate 
real factor at every time step until the total kinetic energy contained in the shell 
matches a predetermined value. The energy dissipation is measured by monitoring 

the multiplying factor. . 

The first inviscid simulation was carried out at 64 3 resolution, and forcing was 
used to fix the energy in the low wavenumber shell k = |k| < 2.5. The simulation 
was initiated from a fully resolved NS forced field at the same resolution and with 
fcmax»7 ~ 2, where r] is the Kolmogorov scale, with a microscale Reynolds number 
Re x « 36 (Jimenez ei a l. 1993). It was not carried to a stationary state due to time 
limitations. As expected, energy accumulated immediately at the high wavenumber 
end of the spectrum, which took a characteristic k 2 shape. Several stages of its 
evolution are shown in Fig. 1. The absence of an appreciable numerical viscosity 
was checked by monitoring the skewness of the velocity derivatives, which rose from 
an initial value of approximately -0.5, characteristic of NS turbulence, to a final one 
of ±10 -3 . 

The inability to obtain a dissipative spectrum from the truncated Euler equations 
themselves shows that an extra condition is needed in the simulation. Since two 
possible similarity solutions exist, two conditions have to be used to decide the 
proportion in which they enter in the final solution. A useful model in which 
to explore what practical conditions might be most appropriate is the “diffusion” 
approximation to turbulent energy transfer proposed by Leith (1967). He obtains 
a differential equation for E(k) as a function of time and wave number 

dE/dt + de/dk = -2vk 2 E. (la) 

where the dissipation e appears as an energy flux in spectral space, and is given by 

e = -Jfc 13 / 2 ^-*- 3 £ 3 / 2 . (lb) 

ok 

Equation (lb) was constructed to have the solutions k 2 and k I in the inviscid 
steady state limit and can be derived from an energy shell model with both forward 
and backward scatter, in the limit of infinitely thin spectral shells (Bell and Nelkin, 
1978). The two similarity solutions are 


e = const., E = (2e/13) 2/3 Ar 5/3 , 

(2) 

e = 0, E ~ k 2 . 

(3) 


Since eq. (1) is nonlinear, the final stationary inviscid spectrum is not simply a 
linear combination of both solutions, but 

E = (2f/13) 2/3 (k~ 5/2 + afc 3 ) 2/ . 


(4) 
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Figure 2. Three dimensional energy spectra for Euler simulations in which 
the total energy at high wavenumbers, k > k c , is kept fixed, in addition to the 

forcing at low wave numbers, (a) k c = fr max - 1. : energy in damped shell, 

E c = 7.5 x 10 2 ; : 10~ 3 . (b) k c = 2* ma */3. line: E c = 0.05; : 

approx. 0.027, see text; : 5 x 10 -3 ; : 5 x 10 -4 


This solution has the right Kolmogorov form up to an (arbitrary) constant but 
contains an extra term coming from the equipartition component. Two conditions 
are needed to adjust the dissipation, f, and the constant a. Note that in the model 
as well as in the actual flow, the equipartition spectrum is conservative. Because 
the behaviors of the two spectra are very different, one of them increasing and the 
other decreasing with k, the two constants are determined almost independently by 
disjoint wave number ranges. Fixing the energy in a low wave number shell essen- 
tially fixes the dissipation because k 2 is small in that range, while fixing the energy 
in a high wave number shell does the same for the equipartition component. This is 
consistent with the common statement that the large scales of the flow fix the dis- 
sipation, but it is true only if the equipartition component, a, can be guaranteed to 
be 0(1). The previous inviscid example shows that this only happens if the energy 
at high wave numbers has been forced to be small. The role of viscosity, from the 
point of view of the large scales, can be seen as that of avoiding the accumulation of 
energy at high wavenumbers and thus of controlling the k 2 component of the spec- 
trum. In general, spectra with an insufficient dissipation mechanism look “almost 
right at low wave numbers but have a large “hook” at high wave numbers, where 
the k component takes over. Small spectral hooks of this kind are common in the 
last few wave numbers of direct numerical simulations. 

2.2 Constrained simulations 

Figure 2a shows the results of two attempts at fixing the energy both in the 
high and in the low spectral shells for a direct Euler simulation. The numerical 
parameters and the forcing scheme in the low wave numbers are exactly the same 
as for Fig. 1, but the code also adjusts independently the last spectral shell (of width 
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AJt = 1) by multiplying all the modes in it by a new damping factor, computed 
every time step, that brings the total energy in the shell to the desired level. This 
proves to be effective only when the energy chosen is high enough but fails if it is 
too small, in which case the rest of the spectrum decouples from the damped shell. 

A little thought shows that too small an energy in the high wave numbers inhibits 
the energy transfer and that the next lower shell cannot dump energy into it fast 
enough. This is no more than saying that the Reynolds stresses, and the turbulent 
energy transfer, can only work if there are small scales to generate them, and that 
the dissipation of a laminar flow is much smaller. The result is an accumulation of 
energy similar to the one in the inviscid case. In fact, the limit in which the energy 
in the last shell is set exactly to zero is numerically equivalent to de-aliasing, and 
it would produce the same result as the unconstrained simulations in Fig. 1. 

Controlling the energy over a thicker shell is more efficient since the lower modes 
have more spectral triads into which the energy can be transferred (or more active 
small scales). The result of fixing, in the same way as before, the total enerp 
contained in all the modes with wave numbers k > k c — 2fc max /3 is shown in Fig. 
2b. It is now possible to control the energy of the high modes almost well enough 
to obtain a k~ 5 ^ 3 extrapolation from the modes that sure completely resolved by the 
simulation. The solid line in Fig. 2b is defined in that way; the energy above k c 
is controlled at each step to coincide with the value it would have if its spectrum 
was a smooth inertial extrapolation of E(k c ). That choice seems to be close to 
optimal from the point of view of obtaining a fc~ 5/3 over most of the inviscid scales. 
However, in this case also, an attempt to make the energy too small in the damped 
shell gross contamination of the spectrum with the equipartition component. 

In all these simulations, the resolved scales are inviscid. In Fig. 1, all the spectral 
modes, except those in the lower forced shell, satisfy the Euler equations, and 
since there is no energy dissipation, the solution tends to equipartition. In the 
simulations in Fig. 2, only those modes with k < k c satisfy the inviscid equations. 
They receive energy from the lower forced shell and lose it to those modes above 
k c , where it is removed every time step by the damping factor. This factor acts as 
an effective viscosity which is zero for all modes below k c , but which is active above 
that threshold and which can be considered as an extreme case of wavenumber 
dependent hyperviscosity. 

Fornberg (1977) used a similar simulation scheme to compute two dimensional 
turbulence. He let the enstrophy in a high shell grow for some time under the effect 
of the cascade and zeroed it periodically. The present scheme is probably a better 
representation of the true cascade dynamics in that it gives the flow a stationary 
energy reservoir with which to interact 

In both schemes, the effect of the damped high-wavenumber shell is to provide 
a spectral bucket into which energy can be dumped naturally by the nonlinear 
interactions. It appears from the results of the simulations that a simple uniform 
damping factor, which adjusts the amplitudes of the small scale modes but does 
not modify their phase relationships, provides an approximation to the small scale 
flow that is good enough to produce a reasonably accurate level of energy transfer 
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into the damped shell. The key seems to be having the right amount of energy in 
that shell to produce the correct cascade rate. This suggests that an ideal transfer 
could be obtained if the energy spectrum in the damped shell was adjusted to follow 
exactly the k~ 5 / 3 law. 

A similar experiment was performed by (She &: Jackson 1993, Zhou 1993). Both 
investigators constrained the entire spectrum to follow exactly the k ~ 5 / 3 law and 
studied different aspects of the energy transfer. In both cases, the correct scaling 
behaviors were observed find the Kolmogorov constant fell in the range of the ex- 
perimental values. Shtilman and Chasnov (1992) also studied a constrained Euler 
system in which the spectrum was forced to remain exactly equal to that of a direct 
simulation of the same flow, run in parallel at higher resolution. Their goal was to 
determine whether the statistics of the constrained flow would be similar to those 
of the DNS field at the same time when filtered down to the same resolution. The 
correspondence was shown to be only approximate, especially for the higher statis- 
tical moments, but the low order statistics, those involved in the energy transfer 
mechanism, agreed well. 

Our goal is different. We axe not only trying to understand the degree to which the 
constrained system approximates locally the energy transfer, but also to determine 
whether an unconstrained system can be made to simulate the inviscid inertial 
spectrum by providing it with a downstream spectral “boundary condition” which 
absorbs energy in the right way. If this turns out to be possible, the unconstrained 
system might better approximate the statistics of the large scales, even though 
those of the small scales are only accurate to the extent of providing a correct 
energy sink. Such a scheme also offers a better chance of being adaptable to the 
study of inhomogeneous turbulent flows. We are therefore interested in constraining 
only a partial spectral shell at high wave numbers to a spectrum that is the fc -5 / 3 
extrapolation of the energy in the last fully resolved shell fc c . 

This was done for the dotted line in Fig. 3. The modes in k > k c — &max/2 
were divided into sub- shells of thickness Ak = 1, and an independent damping 
factor was computed for each sub-shell at each time step to bring the spectrum to 
the form E(k , f) = E(k c , t)(k c /k) 5 / 3 . As expected from the previous discussion, 
the behavior of the whole spectrum and not just that of the forced part was much 
closer to inertial than before. When the compensated spectrum is plotted in linear 
coordinates, however, it is seen that the result is still not perfect (Fig. 3b). A perfect 
inertial spectrum in this representation would be strictly equal to the Kolmogorov 
constant. 


2.3 Turbulent viscosity 

What are still missing are the nonlinear interactions involving scales smaller than 
those at the high end of the damped shell, which are not resolved in the simulation. 
Those interactions, linking the unconstrained modes with those in the sub-grid 
scales, act across a spectral gap imposed by the damped shell, and it is generally 
accepted (Kraichnan 1976, Yoshizawa 1982) that under those circumstances the 
effect of the small scales is equivalent to that of a constant eddy viscosity. The 
presence in filtered DNS fields of an eddy viscosity plateau at wave numbers one or 
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FIGURE 3. Three dimensional energy spectra for Euler simulations in which the 
energy spectrum in the high wavenumber shell, k > k mSLX /2 is adjusted continuously 
to fc“ 5 / 3 , in addition to a constant eddy viscosity and to the forcing at low wave 

numbers, (a) Uncompensated, (b) Normalized. : normalized eddy viscosity, 

ot/ao = 0; : 0.8; : 1.4. 


more octaves lower than the filter width has been confirmed experimentally (Lesieur 
and Rogallo 1989). From dimensional considerations (Kraichnan 1976), the viscosity 
coefficient can be written as 


v e =a [E(k m&x )/k m!ix ] 1/2 . (5) 

The normalized eddy viscosity a can be evaluated within the EDQNM approxi- 
mation for infinite Reynolds number as <Xq « 0.267 (Lesieur, 1990, pg. 324). This 
value assumes that the energy spectrum is strictly inertial beyond fc max > with a Kol- 
mogorov constant of 1.4. The effect of adding this viscosity to the previous simula- 
tion is shown in Fig. 3. All modes are first marched according to the NS equations, 
with the viscosity given by rescaled to fc -5 / 3 . Different values of a produce slightly 
different results, but the spectrum seems to fit best a single fc -5 / 3 law for a value 
slightly smaller than the EDQNM result, a as 0.214 = 0.8ao. The corresponding 
Kolmogorov constant is in the range 1.8-1. 9 (see Fig. 3b), in excellent agreement 
with that obtained in (She & Jackson 1993). This value is somewhat larger than the 
usually accepted one of 1.5, and the reason for the discrepancy is not clear. Zhou 
(1993) obtains a value of approximately 1.5 using a method which is conceptually 
similar to the one in (She & Jackson 1993). It has been suggested (Rogallo, pri- 
vate communication) that the discrepancy between the numerical and experimental 
values may result from the relatively small number of numerical modes in the sim- 
ulation (64 3 ) which robs the cascade of the interactions with the low wavenumber 
modes at scales larger than the size of the computational box. If that were so, the 
correct value would be approached for larger simulations. A single experiment in a 
128 3 box, using the same subgrid mode, gave a slightly lower constant in the range 
1.75-1.85, but the difference was too small for a definite conclusion. 




Energy transfer in isotropic turbulence 


179 



Figure 4. Equilibrium three-dimensional energy spectrum for 64 3 numerical 
boxes in which energy is forced both at the low wavenumber shell, and for modes 
above k c = & raax /2. : Total energy in the high shell extrapolated to it -5 / 3 . 

: Energy spectrum in high shell extrapolated. Vertical line is limit of damped 
shell. Turbulent viscosity in both cases: a/a 0 = 0.8. 


In any case, these experiments suggest that, besides the constant eddy viscosity 
introduced to model the interaction with the very small scales, all that is needed to 
correctly model the dissipation at high Reynolds numbers is to have a “sacrificial” 
spectral shell with the right amount of energy, just beyond the scales of interest. 
The nonlinear interactions with the modes in that shell mimic those with scales just 
below the fully resolved ones in the real flow. It is important to realize that there 
is little reason to expect the modes in the damped shell to be correctly computed 

and that they should be filtered out when flow statistics higher than the energy are 
calculated. 

For the infinite Reynolds number case studied here, the optimum damping strat- 
egy seems to be to constrain the spectrum in the sacrificial shell to the extrapolated 
k law. Other similar strategies can probably be developed for flows at lower 
Reynolds numbers. While this is trivial to implement in the present isotropic spec- 
tral code, it is interesting to inquire whether the energy transfer outside the damped 
shell is sensitive to the detailed form of the constrained spectrum, or whether a sim- 
pler strategy might be used. 

Figure 4 shows the results of an experiment in which the total energy in the last 
shell (k > k c = fc max /2), rather than its spectrum, was used as a constraint, with 
the extra addition of the eddy viscosity of Fig. 3. The total energy is adjusted to the 
one that would be contained in the damped layer if its spectrum were extrapolated 
from E(k c ) according to the k 5 / 3 law. The optimum simulation from Fig. 3 is also 

included for comparison, and the results are seen to differ little outside the damped 
layer. 

From a practical point of view, this scheme is much simpler that the full k ~ 5 / 3 
one since it involves only the calculation of a single global energy and of a single 
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damping factor. It could, in principle, be implemented in a finite difference code by 
the use of appropriate filters. 

3. Backscatter 

The results in the previous section are consistent with the classical idea that the 
interaction of a given range of turbulent scales with features of much smaller size 
can be approximately modeled by an eddy viscosity, even though the interaction 
between similar scales can not. A common scheme for the modeling of these near 
interactions is to separate them into an enhanced wave-number dependent eddy 
viscosity and a backscatter component which moves energy from the small to the 
large scales (Kraichnan 1976, Leslie k Quarini 1979). The backscatter component 
is sometimes represented as a random stirring force (Chasnov 1991) because it can 
be shown that the effect of such a force is always to feed energy into the flow. The 
separation into forward and backward cascade has been re-encountered recently in 
a different context. When a direct simulation of turbulence is divided into large and 
small scale components by filtering the velocity field, the flow of energy between the 
two components is found to be randomly positive and negative at different spatial 
locations in the flow, so that the global energy drain to the small scales is only 
the small residue from two large terms of opposite sign (Piomelli et al. 1991). It 
has also been shown that most of this random energy exchange is dominated by 
interactions among structures separated at most by a factor of two in wavenumber 
magnitude (Domaradzki, Liu k Brachet 1993). 

There is no clear correspondence between the two different definitions of backscat- 
ter, and in fact, it can be argued that the existence of regions of positive and negative 
energy transfer has little to do with the presence of a cascade. It is easy to find 
examples of steady laminar flows in which there is no net transfer of energy, but 
which contain regions of localized forward and backward “cascade”. An examina- 
tion of such flows might be useful in understanding more complicated situations. 
We will discuss a simple example which admits an analytical solution. 

Consider a velocity field u which is decomposed into large and small scale com- 
ponents by means of a circular box filter 

G(x) = 1/7T a 2 if |x| < a, G(x) = 0 otherwise. (6) 

Denote filtered, large scale, variables by Em overbar. The subgrid Reynolds stress is 
defined as 

T{j = UiUj — U t Uj. (7) 

It is also possible to define a large scale rate of strain tensor 

Sij = {Ui,j + Uj,i)/ 2 , (8) 

and a local subgrid “Reynolds” energy transfer rate 


Q — TijSij. 


(9) 
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This term appears on the right hand side of the energy equation for the large scales 
and, with the opposite sign, on that for the small scales. Consider now a simple 
steady two dimensional flow formed by a point vortex, of circulation 7, located at 
the origin in the pure straining field 


Ui - -<JX 2 , U 2 = <JX\ . (10) 

It is easy to see that the effect of the filter on the point vortex is to spread it into a 
circular patch of uniform vorticity u> a = y/ 4 ira 2 with radius 2a, while the straining 
flow (10) is left invariant. The full Reynolds stress and rate of strain tensors can be 
computed explicitly (see Appendix) and the dominant terms near the origin are 

r _ 7 2 log(q/r) / sin 2 <j> — sin 4 > cos <f> \ 

13 27 r 2 a 2 \ — sin 0 cos cos 2 <f> y’ 

which is due to the vortex, and 

= T)- < 12 ) 

which is due to the stagnation flow, where (r, <j>) are polar coordinates based on 
(xi, X2). When r <C a, the transfer term behaves as 

7 2 <j 

Q = 2^7 lo g( Q / r ) sin 2 <t>, (13) 

which is not zero pointwise even though its integral vanishes identically when taken 
over the whole plane, as it must for any stationary flow. 

The interpretation is easy. The small scales of the flow lie near the vortex, where 
gradients are large. The vortex alone produces no transfer since fluid particles follow 
circular trajectories at a constant distance from the center and, therefore, remain 
at the same scale . When the stagnation flow is added, the particles move towards 
the vortex when they pass near <f> = n /4 or 5 tt/ 4, and away from it when they pass 
near <j> = 3 n /4 or 77r/4. In the first case, the particle moves from a large scale to a 
small scale environment, and the energy appears to cascade to higher wave numbers. 
In the second case, the opposite is true, and the cascade is reversed. Note that in 
reality there is only a reversible deformation of fluid particles as they are carried 
by the flow into regions in which the dominant local gradients are larger or smaller. 
It can be shown that in situations less symmetric than the one discussed here, the 
energy exchange is most active when the scales of the advecting flow are roughly 
similar to those of the flow generating the gradients. The “true” cascade, resulting 
in an irreversible transfer of energy from one scale to another, comes from those 
interactions in which the two scales are very different and in which the deformations 
are not reversible. The idea that most of the interaction between neighboring scales 
is not a simple energy diffusion process was put forward forcibly by Kraichnan 
(1976), Rose (1977), and many others afterwards. 


182 


J. Jimenez 


The situation is complicated further because the quantity Q is only part of the real 
subgrid energy flux e(x y t) = -D(u 2 /2)/Dt, although it is often used to represent 
the entire flux because the volume integrals of the two are the same. The local 
difference, which contains the divergence of pressure-velocity and stress-velocity 
products, can be of the same order of magnitude as the flux itself, and in fact, it is 
easily seen that the sign of (13) is not even necessarily the same as that of the real 
transfer because it is compensated by the divergence terms. A consequence is that 
while it makes sense to speak of energy transfer over large volumes of the fluid or of 
transfer between widely separated scales in which the surface terms tend to cancel 
statistically, the local energy transfer due to the type of eddy advection displayed 
by the previous example is a poorly defined quantity that may be difficult to model 
in any locally deterministic way. The procedure outlined in the previous section 
avoids this by providing a model for the subgrid scales rather than trying to model 
the stresses that they produce. 

4. Discussion and future work 

We have shown that, at least for the calculation of the energy spectrum, isotropic 
turbulence at high Reynolds numbers can be modeled rather simply by damping 
a high wavenumber spectral shell in such a way that the energy it contains is the 
same as that which would be contained in a corresponding shell within an inertial 
range. The simulations are improved somewhat by the addition of a constant eddy 
viscosity to account for the nonlinear interactions with scales smaller than those 
contained in the constrained shell. The magnitude of the viscosity coefficient is a 
function of the energy in the last spectral shell and agrees approximately with the 
predictions of EDQNM calculations. The interactions with the flow scales in the 
constrained shell are approximately accounted for by providing that shell with the 
right energy, and they include all the backscatter from small to large scales since 
the effect of the eddy viscosity is purely dissipative. 

The effect on the resolved scales of varying the normalized eddy viscosity, the only 
adjustable parameter in the model, is slight (Fig. 3) although it has a large effect on 
the behavior of the damped shell. This is best explained by the fact that most of the 
energy transferred from the resolved scales goes into the spectral octave immediately 
below it, but that the scales themselves are insensitive to the eventual fate of the 
energy once it has been transferred. For example, in both of the simulations in 
Fig. 4, approximately 72% of the energy is lost to the damped shell while only the 
remaining 28% is withdrawn directly from the resolved scales by the eddy viscosity. 
These values are broadly comparable to those predicted by the simplified analysis 
of Tennekes L Lumley (1972), who on the basis of qualitative arguments about the 
physical mechanism of the cascade, conclude that the energy transfer between shells 
with wave numbers k and k f > k should vary as (fc/fc / ) 7 / 3 . This behavior has been 
verified in (Zhou 1993) for shells with large wave number disparities. An integration 
of that law results in a prediction that 61% of the transfer should go to the next 
wave number octave. 

The insensitivity of the resolved scales to the details of the energy dissipation in 
the damped shell can be seen from the comparatively small differences between their 
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spectra in the a = 0 and a/a 0 = 0.8 simulations in Fig. 3. The energy budgets 
of the damped shell in both cases are very different. For a = 0, energy is removed 
only from the damped shell by the damping factor, while for a/a o = 0.8, 28% is 
removed from the resolved scales by viscosity, 52% is removed by viscosity from the 
damped shell, and only the remaining 20% is removed by the damping factor. It is 
this locality of the energy transfer, coupled with insensitivity to the details of the 
unresolved energy shells, that makes modeling practical. 

The present set of experiments belongs to the recently developed class in which 
part of the flow field is computed as auxiliary to the calculation of the largest 
scales and eventually discarded. The best known of them is the Dynamic model in 
(Germano et al. 1991). In that model, the subgrid dissipation is modeled by an 
eddy viscosity whose magnitude is computed by enforcing similarity of two different 
scales of the flow. Because the similarity is applied locally, the presence of pointwise 
backscatter results in locally negative viscosity coefficients that cause numerical 
problems. It follows from the arguments in §3 that such diffusive models for the 
backscatter are unlikely to be successful. 

We have followed a different route in that we have simulated the near interactions 
by providing an actual energy shell of roughly the right amplitude with which 
the flow can interchange energy naturally. There is no guarantee that the flow 
structures whose scales are in this shell are correctly computed, and they should 
be filtered before the results are used. This is true of the smallest scales in most 
LES calculations. On the other hand, the present scheme is intrinsically dissipative 
and unlikely to develop numerical instabilities. As in the Dynamic model, most of 
the dissipation is computed from the flow parameters themselves, and the residual 
eddy viscosity can be computed from analytic arguments. Furthermore, since the 
eddy viscosity scales with the turbulent energy, the laminar limit is contained in 
the model. Although the discussion in this paper is framed in terms of a spectral 
approach, it is indicated at the end of §2 how the model could be adapted to finite 
difference codes. 

The experience gained in this work provides some insight into the reasons for 
the success of the Dynamic model in many situations. The similarity assumption 
in the Dynamic model essentially provides the right eddy viscosity for a smooth 
extrapolation of the spectrum of the resolved scales (the “test” filter) into the 
smaller scales and, therefore, provides a high wavenumber component with the 
right extrapolated energy. We have shown here that this seems to be enough to 
provide the field with the right dissipation. 

The present model has been developed as a research tool for the study of energy 
transfer in isotropic turbulence and has only been tested on low order statistics. 
Since most high order statistical moments are known to be dominated bv the small 
scales of the flow, they are unlikely to be captured correctly by this, or by any, LES 
model. More important, we have only dealt with homogeneous turbulence, and as 
outlined in the introduction, there are reasons to believe that some of the conclusions 
might not be directly applicable to inhomogeneous situations. In particular, the 
crucial independence of the large scale flow from the small scale details may not be 
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true in any local sense for those cases. Clearly, much more work is needed in all 
those directions. 
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Appendix: Energy transfer in a simple stationary flow 

We present here the detailed analysis of the flow discussed in §3, formed by a 
point vortex at the stagnation point of a straining flow in eq. (10). We define 
polar coordinates (r, <f>), but express all vectors and tensors as components on the 
Cartesian axes (a?i, # 2 ). The velocity field due to the vortex is given by 

7 

Uv = 27tr ( sin<?!s “ COS< £)> (A - 1) 

and that due to the straining flow is 


u s = — <rr (sin<£, cos<£). 


(A -2) 



186 


J . Jimenez 


After filtering, the vortex becomes a uniform circular vortex patch. We will only 
be interested in the behavior of the energy inside the patch, where r < 2 a. The 
filtered velocity is 

"yr 

u„ = g-j(sin^,-W), (A -3) 


while the filtered straining flow is equal to itself because of linearity. The Reynolds 
stresses are computed according to equation (7) and contain terms coming from the 
interaction of the vortex with itself, from the interaction of the straining flow with 
itself, and from the mutual interaction, 


T"vv — 


4tt 2 o 2 


Tvs = 


7 a 


r2 1 1 

( sin 2 <t> — sin <j> 

16a 2 J 1 

^ — sin <t> cos <f> cos 2 

[ r 2 1 

-j 

/ — sin 2 <f> 0 \ 

4a 2 . 

y 0 cos 2 <t> ) ’ 


')■ 


2 \-si 


sin 2 <t> 
sin <j> cos <j> 


— sin <f> cos <j> 
cos 2 <£ 


)• 


(A -4) 

(A -5) 
(A -6) 


The dominant term when r <C a is the logarithmic term in (A-4), which is reflected 
in eq. (11) in the body of the paper. In the same way, the filtered rate of strain 
tensor has a term coming from the vortex and another coming from the straining 
field, 

-5 _ ~7 f-sin2<t> cos2 <f>\ ^a_ 7 \ 

“ 167ra 2 V cos2< t> sin 2^ J ’ 1 ’ 


S, = -a 



(A -8) 


Both are of the same order near the origin, but (A-7) is orthogonal to the dominant 
stress, and the vortex does not generate energy transfer when interacting with itself. 
The “Reynolds” transfer, Q, near the origin is given by eq. (13). 

Most of this transfer, however, is spurious and cancels with various divergence 
terms in the energy equation. The real transfer of kinetic energy from the subgrid 
scales into the filtered flow field is 


Du 2 / 2 
f ~ Di 


= u • Vu, 


(A -9) 


and can be computed exactly as 




\647t 2 




r 2 sin 2$. 


(A -10) 


Its angular structure is the same as that of Q, but its sign depends on the balance 
between the increase of the azimuthal velocity as a particle approaches the vortex 
and the decrease of its radial velocity as it approaches the stagnation point. 


